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Applications  to  ocean  acoustic  data  from  a  towed  array  and  to  speech  processing  are 
presented  for  an  improved  optimal  time-domain  beamformer,  which  involves  optimizing  over 
all  possible  source  bearings  and  time  series  for  multiple  sources  using  simulated 
annealing.  The  convergence  of  the  parameter  search  is  accelerated  by  accepting  time  series 
perturbations  only  when  the  energy  decreases.  A  comparison  with  the  conventional 
delay-and-sum  beamformer  illustrates  that  the  optimal  beamformer  handles  larger  receiver 
spacing  and  larger  source-to-receiver  ratio.  Periodic  ambiguities  are  essentially 
eliminated  by  using  irregular  receiver  spacing  and  the  improved  search  algorithm.  Weak 
sources  are  handled  with  fractional  beamforming.  Noise  cancellation  is  possible  if  the 
parameters  of  the  noise  are  included  in  the  search  space.  Two-dimensional  localization 
is  performed  for  nearby  sources. 

PACS  numbers:  43.30.Yj,  43.60.Gk 


INTRODUCTION 

Conventional  beamforming  techniques,  such  as  the 
delay-and-sum  beamformer,  involve  an  ambiguity  function 
that  depends  on  a  single  parameter  corresponding  to 
source  bearing.  The  estimates  of  the  source  bearings  cor¬ 
respond  to  the  peaks  in  the  ambiguity  function.  Consider¬ 
ing  that  conventional  beamforming  techniques  involve  col¬ 
lapsing  a  large  number  of  parameters  (the  source  bearings 
and  time  series)  into  a  single  parameter  (the  steering  pa¬ 
rameter  of  the  ambiguity  function),  it  is  rather  amazing 
that  they  often  perform  well  for  problems  involving  mul¬ 
tiple  sources. 

The  optimal  beamformer  estimates  the  source  bearings 
and  time  series  all  at  once  by  optimizing  an  energy  function 
that  depends  on  all  of  the  source  bearings  and  time  series. 1 
This  beamformer,  which  is  a  generalization  of  a  frequency- 
domain  beamformer,2  is  practical  with  simulated  anneal¬ 
ing,  an  efficient  Monte  Carlo  method  for  optimization 
problems  involving  large  numbers  of  parameters.3,4  By 
working  with  all  of  the  unknowns,  the  optimal  beamformer 
easily  utilizes  and  benefits  from  a  priori  information. 1  The 
only  approach  for  utilizing  a  priori  information  with  con¬ 
ventional  beamforming  is  to  define  a  new  ambiguity  func¬ 
tion,  which  by  conservation  is  likely  to  inhibit  performance 
in  some  way.  In  this  paper,  we  improve  the  search  algo¬ 
rithm  for  the  optimal  beamformer,  demonstrate  that  the 
optimal  beamformer  permits  a  larger  source-to-receiver  ra¬ 
tio  than  the  delay-and-sum  beamformer,  and  apply  the  op¬ 
timal  beamformer  to  data. 

To  implement  the  optimal  beamformer  numerically, 
the  source  time  series  are  discretized.  This  often  amounts 
to  thousands  of  unknown  parameters.  Since  the  energy 
function  passes  through  a  unique  minimum  as  one  of  the 
discretized  time  series  points  is  varied,  an  improved  simu¬ 


lated  annealing  algorithm  is  obtained  by  accepting  time 
series  perturbations  only  if  the  energy  is  decreased.  Since 
the  energy  function  may  pass  through  local  minima  as  one 
of  the  bearings  is  varied,  the  usual  acceptance  criterion  of 
simulated  annealing  is  used  for  the  bearings:  A  perturba¬ 
tion  is  always  accepted  if  the  energy  decreases  and,  to  allow 
escape  from  local  minima,  is  accepted  according  to  the 
Boltzmann  probability  distribution  if  the  energy  increases. 
An  example  is  presented  in  Sec.  II  that  illustrates  the  ac¬ 
celerated  convergence  of  the  improved  algorithm. 

In  Sec.  Ill,  the  optimal  beamformer  is  compared  with 
the  delay-and-sum  beamformer.  Although  it  is  not  possible 
in  general  to  place  a  limit  on  the  number  of  sources  that 
can  be  handled  by  a  particular  array  (e.g.,  bearings  may  be 
determined  for  an  unlimited  number  of  cw  sources  of  dif¬ 
ferent  frequencies),  the  delay-and-sum  beamformer  typi¬ 
cally  requires  several  times  the  number  of  receivers  that  the 
optimal  beamformer  requires.  An  example  is  presented  to 
illustrate  that,  for  a  given  array  of  receivers,  the  optimal 
beamformer  may  perform  well,  while  the  delay-and-sum 
beamformer  completely  fails  with  large  false  peaks.  In  Sec. 
IV,  we  illustrate  an  advantage  of  using  irregularly  spaced 
receivers.  Combined  with  the  new  search  algorithm,  irreg¬ 
ular  spacing  essentially  eliminates  time  series  ambiguities 
associated  with  periodic  functions. 

A  method  that  we  refer  to  as  fractional  beamforming  is 
described  in  Sec.  V.  With  this  approach,  the  optimal  beam- 
former  is  first  applied  to  estimate  the  bearings  and  time 
series  of  the  most  intense  sources.  These  signals  are  sub¬ 
tracted  from  the  data,  and  the  reduced  data  are  then 
searched  for  weaker  sources.  If  the  intense  signals  domi¬ 
nate  the  received  time  series,  this  approach  can  be  more 
effective  than  searching  for  all  of  the  sources  at  once.  In 
Sec.  VI,  we  generalize  the  search  algorithm  to  cancel  noise, 
which  is  possible  if  something  is  known  about  the  nature  of 
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FIG.  1.  An  example  configuration  of  the  array  and  sources:  Four  plane- 
wave  time  series  arriving  at  a  linear  array  of  five  hydrophones  from  dif¬ 
ferent  bearings. 

the  noise.  In  Sec.  VII,  we  generalize  the  optimal  beam- 
former  to  the  case  of  incident  spherical  waves  for  two- 
dimensional  localization. 

The  performance  of  the  optimal  beamformer  is  robust 
for  data.  A  speech  processing  application  of  optimal  beam¬ 
forming  is  presented  in  Sec.  VIII.  This  problem  involves  a 
very  large  number  of  parameters  and  complex  signals.  The 
voices  of  speakers  in  a  noisy  crowd  are  isolated,  and  the 
recovered  voices  sound  free  of  interference  from  the  other 
voices.  In  Sec.  IX,  the  optimal  beamformer  is  applied  to 
data  from  an  array  towed  in  the  Atlantic  Ocean.  The  ap¬ 
proach  of  Sec.  VI  is  used  to  cancel  noise  from  this  data. 
The  optimal  beamformer  replicas  are  surgically  removed 
from  the  data,  and  a  conventional  beamformer  is  used  to 
process  the  reduced  data. 


I.  THE  OPTIMAL  BEAMFORMER 

As  illustrated  in  Fig.  1,  a  linear  array  of  sensors  re¬ 
ceives  the  M  plane- wave  acoustic  signals  pm(t)  from  the 
directions  0m.  A  subset  of  the  source  bearings  and  each  of 
the  signal  time  series  are  unknown.  The  signal  Pn{t)  re¬ 
ceived  by  the  nth  receiver  is 

M 

PnV)=  £  sin  0m) +¥„(*).  (1) 

m=  l 

The  delays  rn=xn/c  are  defined  in  terms  of  the  hydro¬ 
phone  locations  xn  and  the  sound  speed  c.  The  noise  term 
¥„(*)  may  be  due  to  additional  sources,  ambient  noise,  and 
other  types  of  noise.  The  replica  signal  Qn{t)  received  by 
the  «th  receiver  is  defined  by 


M 

QnU)=  £  qmV+Tn  sin  <f>m),  (2) 

m=  1 

where  the  test  parameters  qm(t)  and  <f>m  are  the  test  time 
series  and  test  bearings.  The  optimal  beamformer  estimate 
for  pm(t)  and  Qm  is  obtained  by  minimizing  the  energy, 

N  r 

=  £  I  [^n(0  —  £?«(*) ]2  dt>  (3) 

n=  1  J 

over  the  test  parameters. 

Each  of  the  source  time  series  is  discretized  so  that  the 
integral  in  Eq.  (3)  is  approximated  by  a  sum.  Hundreds  or 
thousands  of  points  are  typically  used  to  represent  each 
time  series.  The  simulated  annealing  algorithm  involves 
perturbing  each  parameter  one  at  a  time  and  evaluating 
A E,  the  change  in  E.  When  one  of  the  time  series  points  is 
perturbed,  only  a  few  terms  in  the  discretized  integral  are 
affected.  Since  A E  can  be  computed  efficiently,  optimal 
beamforming  is  an  excellent  application  of  simulated  an¬ 
nealing.  Although  all  of  the  terms  in  the  sum  are  affected 
when  a  source  bearing  is  perturbed,  this  does  not  compro¬ 
mise  efficiency  because  the  number  of  source  bearings  is 
much  smaller  than  the  number  of  time  series  points. 

The  solution  that  minimizes  E  may  not  be  unique.  For 
the  case  of  evenly  spaced  receivers  x„=n  Ax  and  two 
sources,  for  example,  E  vanishes  for 


0i(O=/>i(O +/(*)» 

(4) 

** 

< 

*** 

II 

<N 

(5) 

where  /(t)  is  an  arbitrary  periodic  function  of  period  (Ajc/ 
c)(sin  d2— sin  0j).  This  ambiguity,  which  involves  multi¬ 
ple  terms  for  problems  involving  more  than  two  sources, 
was  suppressed  in  Ref.  1  by  assuming  that  the  time  series 
have  compact  support.  This  assumption  is  usually  not  valid 
for  applications. 

During  the  parameter  search,  the  ambiguity  does  not 
present  itself  while  the  bearing  parameters  are  wandering. 
After  two  or  more  of  the  bearings  lock  in,  however,  the 
ambiguous  function(s)  may  grow  without  bound.  In  terms 
of  the  parameter  landscape,  this  is  analogous  to  entering  a 
flat  valley  surrounded  by  mountains  that  get  higher  with 
distance  into  the  valley.  This  is  a  serious  problem  because 
it  is  very  difficult  for  simulated  annealing  to  find  a  way  out 
of  this  type  of  multidimensional  valley.  Fortunately,  it  is 
very  unlikely  to  have  a  first  encounter  with  such  a  valley  at 
a  point  deep  within  because  the  signal  time  series  tend  to  be 
uncorrelated  when  the  bearings  are  wandering. 

We  illustrate  in  Sec.  Ill  that,  in  principle,  all  but  the 
harmless  dc  component  of  the  periodic  ambiguity  may  be 
eliminated  by  using  irregular  receiver  spacing.  Since  the 
time  series  are  discretized,  however,  the  high-frequency 
components  are  not  entirely  eliminated  when  using  the 
simulated  annealing  algorithm  of  Ref.  1.  When  using  the 
improved  search  algorithm  described  in  Sec.  II,  the  high- 
frequency  components  are  also  essentially  eliminated. 
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II.  AN  IMPROVED  OPTIMIZATION  ALGORITHM 


Simulated  annealing  is  an  efficient  method  for  solving 
optimization  problems  involving  local  minima  and  large 
numbers  of  unknowns.  For  specific  applications,  it  is  usu¬ 
ally  possible  to  tune  this  method  to  improve  performance. 
For  the  beamforming  problem,  significant  improvements 
were  obtained  by  treating  the  bearings  and  time  series  as  a 
mixture  of  two  parameter  types. 1  The  energy  change  due  to 
a  perturbation  of  a  point  in  the  discretized  time  series  is 
usually  much  smaller  than  the  energy  change  due  to  a 
perturbation  in  the  bearing  of  a  source.  If  one  were  to  treat 
all  of  the  parameters  the  same,  the  bearings  would  freeze 
out  of  the  mixture  (usually  at  the  wrong  values)  at  a  rel¬ 
atively  high  temperature.  Using  the  analogy  of  the  anneal¬ 
ing  of  a  crystal,  this  difficulty  was  overcome  by  scaling  out 
the  freezing-point  difference  in  the  acceptance  probabili¬ 
ties. 

In  this  section,  we  present  an  improved  simulated  an¬ 
nealing  algorithm  in  which  perturbations  of  the  time  series 
parameters  are  accepted  only  when  the  energy  decreases. 
This  approach  is  valid  because  the  energy  function  is  a 
parabola,  which  has  a  unique  minimum,  for  each  time  se¬ 
ries  parameter.  As  with  the  simulated  annealing  algorithm 
of  Ref.  1,  each  parameter  is  perturbed  once  during  each 
iteration,  one  parameter  at  a  time.  As  described  in  Ref.  1, 
the  Boltzmann  acceptance  probability  is  used  for  the  bear¬ 
ing  perturbations,  which  are  selected  using  a  cubic  distri¬ 
bution  to  allow  the  fast  simulated  annealing  cooling  sched¬ 
ule.5  The  time  series  perturbations  are  selected  with  the 
approach  described  in  Ref.  1,  but  only  those  perturbations 
that  lower  the  energy  function  are  accepted. 

The  improved  algorithm  suppresses  the  high- 
frequency  components  of  the  periodic  ambiguity,  which  is 
described  in  Sec.  I  for  the  two-source  case.  With  the  orig¬ 
inal  algorithm,  the  ambiguous  function  does  not  begin  to 
grow  until  after  the  bearings  have  locked  in  (not  necessar¬ 
ily  to  the  correct  values)  because  the  ambiguous  signal 
components  cancel  each  other  only  for  a  particular  set  of 
bearings.  With  the  original  algorithm,  the  ambiguous  func¬ 
tion  may  continue  to  grow  after  a  nearly  optimal  state  is 
reached,  and  the  recovered  signals  usually  contain  higher 
frequencies  than  the  data.  With  the  improved  algorithm, 
the  ambiguous  function  does  not  grow  after  the  bearings 
lock  in  and  the  energy  function  reaches  a  minimum,  and 
the  recovered  signals  tend  to  be  as  smooth  as  the  data.  In 
terms  of  the  parameter  landscape,  the  improved  algorithm 
tends  not  to  penetrate  deeply  into  the  flat  valleys  described 
in  Sec.  I  because  it  only  accepts  downhill  perturbations  in 
the  time  series. 

We  illustrate  the  performance  of  the  improved  simu¬ 
lated  annealing  algorithm  with  example  A,  which  was  con¬ 
sidered  in  Ref.  1  in  terms  of  dimensionless  variables.  To 
repeat  this  example  in  dimensioned  variables,  we  take 
c=1000  m/s  and  work  in  meters  and  milliseconds.  Since 
the  signals, 


TABLE  I.  Data  for  the  examples:  c  is  the  speed  of  sound  in  m/s;  m  is  the 
source  index;  6m  is  the  mth  source  bearing;  Am  is  the  /nth  source  ampli¬ 
tude;  and  am  is  the  time  lag  in  ms  for  the  mth  Gaussian  source,  which  is 
described  by  Am  and  wm  in  ms. 
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are  of  compact  support,  we  assume  compact  support  in  the 
beamforming  algorithm.  The  parameter  values  for  the  five 
sources  are  given  in  Table  I.  We  take  Ajc=20  m  for  the 
seven  equally  spaced  receivers.  The  received  data  appear  in 
Fig.  2,  and  results  appear  in  Figs.  3  and  4  for  the  original 
and  improved  algorithms.  The  improved  algorithm  con¬ 
verges  to  the  correct  source  bearings  several  times  faster 
than  the  original  algorithm  for  this  problem.  The  improved 
algorithm  accurately  recovers  the  source  time  series.  Irreg¬ 
ular  receiver  spacing  is  applied  to  this  problem  in  Sec.  IV. 


III.  THE  CONVENTIONAL  BEAMFORMER 


In  this  section,  we  compare  the  performance  of  the 
optimal  beamformer  and  the  delay-and-sum  beamformer, 
which  is  defined  by 


B{6)  = 


J  (  X  ^(f-r,,  sin  0)j 


dt, 


(7) 


where  N^M  and  the  steering  parameter  0  corresponds  to 
source  bearing.  For  the  case  'J/„(r)  =0,  we  substitute  Eq. 
(1)  into  Eq.  (7)  to  obtain 


,  N  M  .  \2 

B{0)=  X  X  pm[t-Tn{sm  0-sin  6m)]  \  dt. 

J  \  n= 1  m= 1  / 

(8) 
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beamformer  with  seven  equally  spaced  receivers  and  five  sources. 


Unless  9^0m  for  some  m,  the  time  series  in  Eq.  (8)  add 
incoherently,  and 


B(d)=  J  [0(Nl/2)]2dt=0(N). 


0  1000  2000  3000  4000 


For  0=6m,  the  mth  time  series  terms  add  coherently,  and 


(b) 


Iterations 


B(6m)=  J  [NpmU)  +  0(N'/2)]2dt=0(N2).  (10) 

Keeping  the  dominant  term  in  Eq.  (10),  we  obtain 


FIG.  3.  Results  for  example  A,  which  involves  the  optimal  beamformer 
with  seven  equally  spaced  receivers  and  five  sources,  (a)  The  original 
simulated  annealing  algorithm  converges  after  approximately  3000  itera¬ 
tions,  (b)  The  improved  simulated  annealing  algorithm  converges  after 
approximately  500  iterations. 


B(ej~N2  J  [pm(t)?dt.  (11) 

For  example  B,  we  apply  the  delay-and-sum  beam- 
former  to  a  problem  involving  the  source  bearings  and  time 
series  of  example  A.  The  estimates  of  B(0m)  given  by  Eq. 
(11)  are  approximately  proportional  to  A2mwm  for  these 
sources.  We  consider  four  equally  spaced  arrays:  array  1 
(N=l,  Ax =20  m);  array  2  (iV=  14,  Ax =20  m);  array  3 
(N=  14,  Ajc=10  m);  and  array  4  (iV=28,  Ax=  10  m). 
Results  for  example  B  appear  in  Fig.  5.  For  array  1,  which 


was  used  successfully  with  the  optimal  beamformer  in  ex¬ 
ample  A,  the  delay-and-sum  beamformer  has  large  false 
peaks  near  broadside  and  endfire.  The  false  peaks  are 
smaller  for  array  2.  Although  there  are  no  major  false 
peaks  for  array  3,  two  of  the  weaker  source  peaks  are 
barely  separated  from  sidelobes.  All  five  of  the  peaks  are 
well  defined  for  array  4  with  levels  approximately  equal  to 
the  values  predicted  by  Eq.  (11).  Example  B  illustrates 
that  the  optimal  beamformer  can  handle  larger  receiver 
spacing  and  larger  source-to-receiver  ratio  than  the  delay- 
and-sum  beamformer. 
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IV.  IRREGULAR  RECEIVER  SPACING 
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FIG.  4.  Results  for  example  A,  which  involves  the  optimal  beamformer 
with  seven  equally  spaced  receivers  and  five  sources.  The  time  series 
recovered  with  the  improved  algorithm  (solid  curves)  agree  with  the  true 
time  series  (dashed  curves). 


The  performance  of  a  beamforming  algorithm  is  sen¬ 
sitive  to  the  array  parameters.  For  example,  the  results  of 
Sec.  Ill  illustrate  the  importance  of  N  and  Ax.  These  are 
not  the  only  parameters  that  can  be  varied.  For  example, 
some  beamforming  methods  can  be  enhanced  significantly 
by  selecting  the  receiver  locations  so  that  the  array  samples 
a  wide  range  of  length  scales.6  In  this  section,  we  illustrate 
the  advantage  of  using  irregular  receiver  spacing  with  the 
optimal  beamformer. 

For  example  C,  we  consider  the  sources  of  example  A 
and  drop  the  assumption  that  the  source  time  series  have 
compact  support.  When  the  bearings  are  assumed  to  be 
known  a  priori ,  an  evenly  spaced  array  with  N= 9  and 
Ax  =  20  m  recovers  the  source  time  series  appearing  in  Fig. 
6,  which  are  corrupted  by  the  periodic  ambiguity.  Results 
appear  in  Fig.  7  for  an  array  of  nine  receivers  with  spacings 
(from  the  end  of  the  array  corresponding  to  positive  bear¬ 
ing)  of  2,  17,  5,  23,  7,  19,  3,  and  13  m.  We  used  two  more 
receivers  than  for  example  A  because  some  of  the  receiver 
spaces  are  small.  With  this  array,  the  source  bearings  are 
recovered  and  the  recovered  time  series  agree  with  the  true 
time  series  with  only  small  hints  of  the  periodic  ambiguity. 


FIG.  5.  Results  for  example  B,  which  involves  the  conventional  delay-and-sum  beamformer  with  equally  spaced  receivers  and  five  sources,  (a)  There 
are  large  false  peaks  near  broadside  and  endfire  for  array  1.  (b)  The  levels  of  the  false  peaks  are  reduced  for  array  2.  (c)  Although  there  are  no  major 
false  peaks  for  array  3,  two  of  the  source  peaks  are  barely  separated  from  sidelobes.  (d)  The  five  source  peaks  are  well  defined  for  array  4. 
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FIG.  6.  Results  for  example  C,  which  involves  the  five  sources  of  example 
A  without  the  compact  support  assumption.  Using  nine  evenly  spaced 
receivers  and  assuming  that  the  source  bearings  are  known  a  priori,  the 
recovered  time  series  are  corrupted  by  the  periodic  ambiguity. 


V.  FRACTIONAL  BEAMFORMING 

All  beamforming  methods  have  difficulty  determining 
the  bearings  of  weak  sources  that  are  dominated  by  inter¬ 
ference  from  strong  sources.  In  this  section,  we  present  an 
effective  approach  for  handling  weak  sources  with  the  op¬ 
timal  beamformer.  We  refer  to  the  approach  as  fractional 
beamforming,  which  involves  two  steps.  The  optimal 
beamformer  is  first  applied  to  estimate  the  bearings  and 
time  series  of  the  strong  sources.  If  a  sufficiently  low  energy 
is  achieved  (this  occurs  only  if  the  weak  sources  truly  are 
weak),  the  strong  sources  are  then  subtracted  from  the 
data.  The  optimal  beamformer  is  applied  to  the  reduced 
data  to  estimate  the  bearings  and  time  series  of  the  weak 
sources.  Fractional  beamforming  often  performs  much  bet¬ 
ter  than  the  approach  of  applying  the  optimal  beamformer 
to  search  for  all  of  the  sources  at  once.  The  performance  is 
best  when  the  weak  sources  are  substantially  below  the 
strong  sources.  Similar  behavior  occurs  with  an  approach 
for  extracting  weak  signals  from  noise.7 

For  example  D,  we  consider  a  problem  involving  two 
strong  sources,  one  weak  source  (20  dB  below  each  of  the 
strong  sources),  and  the  irregularly  spaced  array  of  exam¬ 
ple  C.  The  parameter  values  for  the  sources  are  given  in 
Table  I.  The  results  of  applying  the  optimal  beamformer  to 
estimate  the  bearings  and  time  series  of  the  strong  sources 
from  the  raw  data  appear  in  Fig.  8.  Both  of  the  bearings 
and  time  series  are  recovered  accurately  modulo  a  small  dc 
component.  The  results  of  applying  the  optimal  beam- 
former  to  the  reduced  data  appear  in  Fig.  9.  The  bearing  of 


FIG.  7.  Results  for  example  C,  which  involves  the  five  sources  of  example 
A  without  the  compact  support  assumption.  Using  an  irregularly  spaced 
array  with  nine  receivers,  the  bearings  converge  to  the  correct  values  (a), 
and  the  recovered  time  series  (solid  curves)  and  true  time  series  (dashed) 
curves  are  in  agreement  (b). 


the  weak  source  is  accurately  determined.  The  recovered 
time  series  for  the  weak  source  is  also  fairly  accurate.  The 
results  of  applying  the  optimal  beamformer  to  search  for 
all  three  sources  at  once  appear  in  Fig.  10.  Although  the 
bearings  of  the  strong  sources  are  recovered,  the  bearing 
parameter  for  the  weak  source  wanders  aimlessly. 

VI.  NOISE  CANCELLATION 

The  optimal  beamformer  may  be  generalized  to  cancel 
noise  if  an  appropriate  property  of  the  noise  is  known.  In 


J.  Acoust  Soc.  Am.,  Vol.  93,  No.  4,  Pt.  1,  April  1993 


Collins  et  al.:  Optimal  beamforming  1856 


1856 


0  500  1000  1500  2000 

(a)  Iterations 


FIG.  8.  Results  for  example  D,  which  involves  two  strong  sources  and 
one  weak  source.  The  optimal  beamformer  is  applied  to  the  raw  data 
searching  for  two  sources.  For  the  strong  sources,  (a)  the  bearings  are 
recovered,  and  (b)  the  recovered  time  series  (solid  curves)  agree  with  the 
true  time  series  (dashed  curves)  modulo  a  small  dc  component  of  the 
ambiguity. 


this  section,  we  modify  the  optimal  beamformer  to  cancel 
electromagnetic  noise.  This  type  of  noise  usually  consists  of 
discrete  frequencies.  Although  this  type  of  noise  can  be 
removed  using  the  Fourier  transform,  parts  of  the  signals 
may  also  be  removed  with  this  approach.  The  cancellation 
approach,  which  may  be  applied  to  various  types  of  noise, 
removes  noise  surgically. 

The  noise  term  is  assumed  to  be  of  the  form, 


^«(0=v„sin(fi)f+6),  (12) 

where  co  is  the  circular  frequency  of  the  source.  Unless  the 
noise  on  different  receivers  arises  from  different  processes, 
the  unknown  phase  6  is  independent  of  ft  (modulo  ir) 
because  electromagnetic  signals  propagate  much  faster 
than  acoustic  signals.  The  phase  may  differ  by  tt  on  differ¬ 
ent  receivers  (we  have  observed  this  in  data)  if  (for  exam¬ 
ple)  the  noise  is  generated  by  a  varying  magnetic  field  and 


FIG.  9.  Results  for  example  D,  which  involves  two  strong  sources  and 
one  weak  source.  The  optimal  beamformer  is  applied  to  the  reduced  data 
searching  for  one  source.  For  the  weak  source,  (a)  the  bearing  is  recov¬ 
ered,  and  (b)  the  recovered  time  series  (solid  curves)  agrees  with  the  true 
time  series  (dashed  curves). 
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Iterations 

FIG.  10.  Results  for  example  D,  which  involves  two  strong  sources  and 
one  weak  source.  The  optimal  beamformer  is  applied  to  the  raw  data 
searching  for  three  sources.  The  bearings  are  recovered  for  the  strong 
sources  but  not  for  the  weak  source. 
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FIG.  11.  Data  and  results  for  example  E,  which  involves  cw  noise  from 
electromagnetic  interference.  The  best  replica  (dashed  curve)  agrees  with 
the  data  (solid  curve).  The  highest  noise  level  occurs  on  receivers  1 
and  5. 


the  wires  corresponding  to  different  receivers  are  coiled  in 
opposite  directions. 

The  noise  amplitudes  v„  are  assumed  to  depend  on  n 
because  some  receivers  may  be  noisier  than  others.  To  ac¬ 
count  for  phase  reversals,  we  allow  the  sign  of  vn  to  depend 
on  n.  This  approach  may  be  applied  to  cancel  noise  con¬ 
sisting  of  several  frequencies.  The  solution  is  ambiguous  if 
one  of  the  plane  waves  is  incident  at  broadside.  If  at  least 
one  of  the  receivers  is  noise-free,  however,  this  ambiguity 
may  be  eliminated  by  assuming  the  noise  amplitude  van¬ 
ishes  for  that  receiver. 

We  illustrate  noise  cancellation  with  example  E,  which 
involves  the  nine  irregularly  spaced  receivers  used  in  ex¬ 
ample  C  and  the  three  sources  described  in  Table  I.  The 
50-Hz  noise  amplitudes  on  the  receivers  are  2.0,  — 1.0,  1.4, 
0.2,  2.0,  0.5,  —1.5,  1.0,  and  1.2.  The  noisy  array  data  and 
the  replica  time  series  corresponding  to  the  lowest  energy 
state  encountered  appear  in  Fig.  11.  Other  results  for  ex¬ 
ample  E  appear  in  Fig.  12.  All  three  of  the  source  bearings 


(b)  Time  (ms) 


FIG.  12.  Results  for  example  E,  which  involves  cw  noise  from  electro¬ 
magnetic  interference,  (a)  The  source  bearings  are  recovered,  (b)  The 
source  time  series  are  recovered. 


FIG.  13.  Results  for  example  F,  which  involves  two-dimensional  local¬ 
ization  of  two  point  sources  located  along  the  same  bearing. 
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FIG.  14.  The  phrases  for  example  G,  which  involves  speech  processing 
with  three  speakers  and  seven  receivers.  The  phrases  (from  the  top  down) 
are  “minus  forty,”  “minus  twenty,”  and  “plus  fifteen.” 

are  recovered  accurately.  The  source  time  series  and  the 
noise  parameters  are  also  recovered. 


lengths  away.  For  this  problem,  spherical  spreading  is  in¬ 
cluded  and  the  delays  appearing  in  Eqs.  (1)  and  (2)  are 
modified  to  be 


where  x  is  the  source  location  and  xn  is  the  location  of  the 
nth  receiver.  The  coordinate  origin  is  taken  to  be  the  cen¬ 
troid  of  the  array.  If  the  receivers  are  collinear,  it  is  possi¬ 
ble  to  determine  the  two  source  coordinates  corresponding 
to  range  r  and  bearing  0. 

Example  F  involves  two  sources  located  on  the  same 
bearing  6=  —30  deg  at  r=60  m  and  r=  120  m.  We  use  an 
array  of  seven  equally  spaced  receivers  with  Ax=  10  m  and 
the  source  parameters  appearing  in  Table  I.  Results  for 
example  F  appear  in  Fig.  13.  The  algorithm  converges 
quickly,  and  the  bearings  are  recovered  accurately.  The 
range  of  the  closer  source  is  recovered  with  more  certainty 


VII.  SPHERICAL  WAVES 

In  this  section,  we  modify  the  optimal  beamformer  to 
handle  incident  spherical  waves.  This  approach  is  applica¬ 
ble  when  some  of  the  sources  are  less  than  about  ten  array 


Time  (ms) 


FIG.  16.  Results  for  example  G,  which  involves  speech  processing  with 


three  speakers  and  seven  receivers,  (a)  The  bearings  of  the  speakers, 


FIG.  15.  The  array  data  for  example  G,  which  involves  speech  processing  which  are  the  same  as  the  phrases  (in  degrees),  are  recovered,  (b)  The 


with  three  speakers  and  seven  receivers. 


recovered  phrases  agree  with  the  true  phrases  appearing  in  Fig.  14. 
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FIG.  17.  Results  for  the  first  data  segment  of  example  H,  which  involves 
data  from  the  Atlantic  Ocean,  (a)  The  array  data  (dashed  curve)  and  the 
best  replica  (solid  curve).  Note  the  180-Hz  noise  on  receivers  5  and  6.  (b) 
The  energy  and  the  source  bearings,  (c)  The  recovered  time  series  for  the 
tow  ship  (source  1),  the  passing  ship  (source  2),  and  the  distant  ship 
(source  3). 


because  its  wave-front  curvature  is  easier  to  resolve  with 
this  array.  The  source  range  becomes  a  weak  attractor 
(i.e.,  the  parameter  does  not  tend  to  remain  near  the  true 
value  during  the  search)  as  it  increases.  This  beamformer 
reduces  to  the  plane-wave  optimal  beamformer  in  the  limit 
of  large  source  range. 


VIII.  SPEECH  PROCESSING 

In  this  section,  we  apply  the  optimal  beamformer  to 
example  G,  which  is  a  speech-processing  problem  involv¬ 


ing  complex  signals  and  more  than  70  000  unknowns.  The 
phrases  “minus  forty,”  “minus  twenty,”  and  “plus  fifteen” 
were  digitized  using  a  NeXT  computer.  These  phrases, 
which  appear  in  Fig.  14,  were  each  recorded  separately  for 
1.5  s  using  a  16-kHz  sampling  rate.  The  phrases  were  then 
superimposed  using  Eq.  ( 1 )  to  simulate  data  on  an  evenly 
spaced  array  with  N=1  and  Ax =0.5  m.  The  array  data 
appear  in  Fig.  15.  These  jumbled  phrases  sound  similar 
from  receiver  to  receiver.  The  bearings  used  for  the  speak¬ 
ers  (in  degrees)  correspond  to  the  phrases.  We  assumed 
that  c=300  m/s. 
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Although  the  array  data  are  quasisynthetic,  this  prob¬ 
lem  differs  from  the  simulations  presented  in  previous  sec¬ 
tions  in  that  the  signals  are  highly  complex  and  the  number 
of  parameters  is  much  larger.  The  results  of  the  simulated 
annealing  algorithm  appear  in  Fig.  16.  The  bearings  of  the 
speakers  are  recovered  accurately,  and  the  recovered 
phrases,  which  sound  nearly  identical  to  the  original 
phrases,  are  similar  to  the  true  phrases.  Phrase  1,  which 
contains  more  energy  than  the  other  phrases,  is  recovered 
with  the  greatest  accuracy. 

Due  to  the  large  number  of  time  series  points,  it  took 
several  hours  to  solve  this  problem  on  a  SiliconGraphics 
Iris  workstation.  Since  the  optimal  beamformer  is  parallel 
in  the  time  series  parameters,  however,  run  time  on  a  par¬ 
allel  processing  computer  would  be  essentially  independent 
of  the  number  of  time  series  points. 


IX.  OCEAN  ACOUSTIC  DATA 


In  this  section,  we  apply  the  optimal  beamformer  to 
data  from  a  64-receiver  array  towed  behind  the  USNS 
Lynch  in  a  region  of  the  Atlantic  Ocean  south  of  Nova 
Scotia.  The  array  consists  of  four  nested  subarrays  with 
even  spacings  of  1.25,  2.5,  5,  and  10  m.  The  system  was  not 
designed  for  the  optimal  beamformer.  The  receiver  spac¬ 
ings  were  not  optimized  to  suppress  the  periodic  ambigu¬ 
ity.  Some  receivers  exhibited  a  large  amount  of  180-Hz 
electromagnetic  noise,  which  is  not  a  problem  for  the 
frequency-domain  processors  for  which  the  array  was  de¬ 
signed.  We  took  2-s  data  segments  every  half-hour  using  a 
2-kHz  sampling  rate  and  an  analog  filter  to  remove  energy 
above  500  Hz.  We  were  able  to  obtain  solutions  with  the 
optimal  beamformer  and  verify  them  in  nearly  real  time 
using  radar  and  conventional  frequency-domain  processing 
with  larger  subarrays.  We  assumed  that  c=1500  m/s  for 
the  data  presented  here,  which  was  taken  on  14  August 
1990  near  42.7  deg  N,  62  deg  W.  The  ocean  is  approxi¬ 
mately  1400  m  deep  in  this  area.  The  LYNCH  towed  the 
array  toward  the  southeast  at  3  kn. 

For  example  H,  we  consider  portions  of  three  consec¬ 
utive  data  segments  (a  total  time  of  about  1  h)  that  contain 
signals  from  the  tow  ship,  a  passing  ship  that  was  sweeping 
in  bearing,  and  a  distant  ship  that  remained  near  aft  end- 
fire.  We  worked  with  an  evenly  spaced  subarray  of  seven 
receivers  with  Ax =2.5  m.  Since  the  180-Hz  noise  was  ob¬ 
served  to  be  weak  on  receivers  1,  2,  and  3,  we  assumed  that 
the  noise  vanishes  on  these  receivers  to  suppress  the  broad¬ 
side  ambiguity  described  in  Sec.  VI.  We  searched  for  three 
sources  applying  the  a  priori  knowledge  that  the  tow  ship 
was  near  endfire.  Since  the  array  was  submerged,  we  al¬ 
lowed  the  tow  ship  to  be  up  to  40  deg  from  endfire.  We 
searched  for  only  two  sources  in  the  third  data  segment 
because  the  passing  ship  and  the  distant  ship  were  both 
near  endfire. 

Results  for  the  first  data  segment  appear  in  Fig.  17. 
The  array  data  and  the  best  replica  (i.e.,  the  replica  corre- 


FIG.  18.  Reduced  array  data  for  the  first  data  segment  of  example  H, 
which  involves  data  from  the  Atlantic  Ocean.  The  180-Hz  noise  on  re¬ 
ceivers  5  and  6  has  been  canceled. 

sponding  to  the  lowest  energy  encountered)  are  in  good 
agreement.  The  180-Hz  noise  levels  on  receivers  5  and  6 
are  relatively  high.  The  recovered  source  bearings  corre¬ 
spond  to  the  true  bearings  of  the  passing  ship  and  the 
distant  ship.  The  bearing  recovered  for  the  tow  ship  indi¬ 
cates  that  the  array  dipped  to  about  25  deg  below  endfire. 
The  recovered  signal  for  the  passing  ship  contains  lower 
frequencies  than  the  other  signals.  All  of  the  receivers  ap¬ 
pear  to  be  free  of  180-Hz  noise  in  the  reduced  data,  which 
consists  of  the  data  minus  the  best  180-Hz  noise  replica 
and  appears  in  Fig.  18.  Results  for  the  second  and  third 
data  segments  appear  in  Figs.  19  and  20.  The  data  and 
replica  time  series  are  in  good  agreement  for  both  of  these 
cases.  For  the  second  segment,  the  passing  ship  is  louder 
than  the  other  ships.  For  the  third  segment,  the  combined 
signal  from  the  passing  ship  and  the  distant  ship  is  weaker 
than  the  tow  ship  signal. 

To  illustrate  the  ability  of  the  optimal  beamformer  to 
surgically  remove  signals  from  data,  we  consider  example 
I,  which  involves  the  second  segment  of  data  and  an  evenly 
spaced  subarray  of  nine  receivers  with  Ax =5  m.  These 
array  parameters  were  selected  so  that  conventional  beam- 
forming  would  perform  well  for  lower  frequencies.  A 
frequency-azimuth  (FRAZ)  diagram  (which  displays 
beamformed  energy  versus  frequency  and  azimuth)  of  the 
raw  data  appears  in  Fig.  21(a).  To  construct  this  diagram, 
0.5  s  of  data  was  used  to  estimate  the  spectrum,  and  the 
Bartlett  processor  was  evaluated  for  181  bearings  and  220 
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FIG.  19.  Results  for  the  second  data  segment  of  example  H,  which  in¬ 
volves  data  from  the  Atlantic  Ocean,  (a)  The  array  data  (dashed  curve) 
and  the  best  replica  (solid  curve),  (b)  The  energy  and  the  source  bear¬ 
ings.  (c)  The  recovered  time  series  for  the  tow  ship  (source  1),  the 
passing  ship  (source  2),  and  the  distant  ship  (source  3). 


frequencies.  The  passing  ship  appears  near  0=40  deg  and 
the  other  ships  appear  near  the  opposite  endfire  positions. 
There  are  numerous  sidelobes  for  the  three  sources.  These 
false  peaks  can  be  distinguished  from  true  source  peaks 
because  their  locations  vary  with  frequency.  Although  the 
electromagnetic  noise  is  evident  at  various  multiples  of  60 
Hz  (especially  180  Hz),  it  is  weak  relative  to  the  broad¬ 
band  acoustic  signals,  and  we  do  not  bother  to  cancel  it  in 
this  subarray. 

Using  the  a  priori  information  gained  from  example  H, 
we  initialized  the  temperature  at  a  low  value  and  restricted 


the  source  bearings  to  narrow  windows.  After  200  itera¬ 
tions,  the  algorithm  produced  results  similar  to  the  results 
of  example  H.  A  FRAZ  diagram  appears  in  Fig.  21(b)  for 
the  reduced  data  obtained  by  subtracting  the  best  tow-ship 
replica  from  the  data.  All  evidence  of  the  tow  ship  is  gone 
including  its  sidelobes.  The  features  corresponding  to  the 
other  ships  do  not  appear  to  be  affected.  Appearing  in  Fig. 
21(c)  is  the  FRAZ  diagram  of  the  data  minus  the  best 
replicas  for  the  tow  ship  and  the  passing  ship.  The  peaks 
corresponding  to  the  distant  ship  and  its  sidelobes  remain 
intact.  The  FRAZ  diagram  obtained  by  subtracting  all 
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FIG.  20.  Results  for  the  third  data  segment  of  example  H,  which  involves 
data  from  the  Atlantic  Ocean,  (a)  The  array  data  (dashed  curve)  and  the 
best  replica  (solid  curve),  (b)  The  energy  and  the  source  bearings,  (c) 
The  recovered  time  series  for  the  tow  ship  (source  1 )  and  the  passing  ship 
and  distant  ship  combined  (source  2). 


three  best  replicas  appears  in  Fig.  21(d).  Evidence  of  a 
fourth  signal  near  0=20  deg  appears  in  the  low-frequency 
region  of  the  diagram. 

X.  CONCLUSION 

An  improved  simulated  annealing  algorithm  has  been 
developed  for  the  optimal  beamformer.  By  accepting  time 
series  perturbations  only  if  the  energy  function  decreases, 
convergence  is  accelerated  and  high-frequency  noise  asso¬ 
ciated  with  the  periodic  ambiguity  is  suppressed.  Frac¬ 
tional  beamforming  is  an  effective  approach  for  handling 


weak  sources.  The  optimal  beamformer  handles  larger  re¬ 
ceiver  spacings  and  more  sources  per  receiver  than  conven¬ 
tional  beamformers,  which  collapse  all  of  the  unknowns 
into  a  single  parameter.  By  using  irregular  receiver  spac¬ 
ing,  the  performance  of  the  optimal  beamformer  may  be 
enhanced  and  the  low-frequency  component  of  the  periodic 
ambiguity  may  be  suppressed.  It  is  possible  to  cancel  cer¬ 
tain  types  of  noise  with  the  optimal  beamformer  if  the 
nature  of  the  noise  is  understood.  The  optimal  beamformer 
performs  well  for  towed  array  data  and  for  speech 
processing. 
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Frequency  (Hz)  Frequency  (Hz) 
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FIG.  21.  Results  for  example  I,  which  involves  data  from  the  Atlantic  Ocean.  FRAZ  diagrams  for  (a)  the  raw  data;  (b)  the  raw  data  minus  the  best 
tow  ship  replica;  (c)  the  raw  data  minus  the  best  tow  ship  and  passing  ship  replicas;  and  (d)  the  raw  data  minus  all  three  replicas.  The  colors  span  20 
dB,  with  red  corresponding  to  high  intensity  and  blue  corresponding  to  low  intensity. 
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